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Runoff  Prediction  Uncertainty  for  Ungauged 
Agricultural  Watersheds 

By  David  M.  Goldman/  Miguel  A.  Marino/  and  Arlen  D.  Feldman/ 

Members,  ASCE 

Abstract:  a  physically  based  stochastic  watershed  model  is  used  to  estimate 
runoff  prediction  uncertainty  for  small  agricultural  watersheds  in  Hastings,  Ne¬ 
braska.  The  stochastic  nature  of  the  model  results  from  postulating  a  probabilistic 
model  for  parameter  estimation  and  input  errors.  The  key  factors,  assumed  to  con¬ 
tribute  to  prediction  uncertainty  are  errors  in  estimating  infiltration  parameters  and 
moisture  conditions  prior  to  a  rainfall  event.  The  error  distributions  for  parameter 
estimates  are  inferred  from  soil  survey  information,  and  the  error  distribution  for 
moisture  conditions  from  a  regression  between  antcr^ent  precipitation  indices  and 
measured  soil  moisture.  Comparison  of  model  predicted  and  observed  errors  dem¬ 
onstrates  that  the  model  is  conservative  in  that  it  is  biased  towards  overprediction 
of  errors. 


Introduction 

A  common  approach  to  deriving  flood-flow-frequency  curves  for  an  un¬ 
gauged  watershed  is  to  simulate  a  design  storm  with  an  event-oriented  wa¬ 
tershed  model.  As  with  most  models,  a  significant  problem  with  this  ap¬ 
proach.  is  model  parameter  estimation.  One  approach  to  the  problem  is  to 
estimate  watershed  model  parameters  from  generally  available  information 
such  as  topographic  maps  and  soil  maps.  Probably  the  most  difficult  part  of 
the  problem  is  estimating  loss  rates  from  this  type  of  information. 

Loss  rate  parameter  estimation  for  an  event-oriented  model  consists  of 
determining  both  the  soil  moisture  condition  prior  to  the  rainfall  event  and 
the  soil  infiltration  parameters.  The  particular  focus  of  this  study  will  be  on 
determining  the  uncertainty  in  estimating  initial  moisture  conditions  from 
antecedent  precipitation  in^ces  and  infiltration  parameters  from  commonly 
available  soil  survey  information,  particle-size  distribution,  and  porosity,  as 
proposed  by  Rawls  and  Brakensiek  (1982). 

An  examination  of  the  runoff  prediction  uncertainty  caused  by  estimating 
loss  rate  parameters  from  simple  indices  such  as  particle-size  distribution 
will  be  investigated  using  a  physically  based  stochastic  (PBS)  watershed  model 
(Klemes  1978).  The  PBS  model  will  be  formulated  as  a  typical  event-ori¬ 
ented  distributed  watershed  model,  except  that  the  loss  rate  parameters  will 
be  considered  stochastic  (Goldman  1987).  The  stochastic  nature  of  the  loss 
rate  parameters  will  be  represented  by  the  loss  rate  parameters’  probability 
distributions.  The  parameter  probability  distribution  represents  the  estimation 
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uncertainty  when  depending  on  typically  available  information  such  as  the 
simple  indices  proposed. 

TTie  PBS  model  prediction  uncertainty  is  represented  by  the  derived  prob¬ 
ability  distribution  for  runoff  volume  and  peak  discharge.  The  derived  i»x>b- 
ability  distribution  gives  the  runoBf  exceedance  probabilities  for  any  partic¬ 
ular  runoff  event.  In  other  words,  the  derived  probability  distributions  give 
the  chance  that  an  observed  runoff  volume  or  peak  discharge  will  exceed  a 
given  amount  for  a  given  rainfall  event  and  antecedent  moisture  condition. 
The  usefulness  of  the  model-estimated  prediction  uncertainty,  and  in  turn 
the  estimated  uncertainty  in  parameter  values,  will  be  verified  by  comparing 
predicted  and  observed  exceedance  probabilities  for  rainfall-runoff  events 
recorded  for  small  agricultural  watersheds  in  Hastings,  Nebraska. 

Test  Watershed  Description 

The  U.S.  Department  of  Agriculture  (USDA)  maintained  wat»isheds  in 
Hastings,  Nebraska  (The  Central  no  date),  where  rainfall-runoff,  soil-mois¬ 
ture,  and  land-use  data  were  collected  from  the  period  1939-1967  and  main¬ 
tained  in  the  REPHLEX  data  base  (REPHLEX  1983).  The  seven  watersheds 
used  in  this  research  were  all  about  4.0  acres  in  size  with  overland  flow 
slopes  ranging  from  2.0  to  7.0%.  The  watershed  soils  are  dominated  by  the 
Hastings  silt  loam  and  silt  clay  loam  horizons.  Although  data  for  a  number 
of  different  management  practices  existed,  only  runoff  events  from  land  sur¬ 
faces  with  straight  row  contouring  were  considered,  because  this  situation 
corresponds  best  to  the  assumptions  that  will  be  made  in  developing  the  PBS 
model. 

There  was  always  a  rain  gage  within  200  ft  of  the  border  of  any  of  die 
watersheds  tested.  Unfortunately,  however,  the  rain  gages  were  not  shielded. 
The  soil  moisture  data  were  obtained  by  analyzing  fist-sized  soil  samples 
taken  from  field  trenches.  The  major  advantage  of  using  die  Hastings  data 
is  that  they  had  a  relatively  long  record,  which  was  important  for  verifying 
the  predictions  of  the  PBS  model. 

Model  Formulation 

A  single-event  approach  to  modeling  runoff  was  chosen  to  avoid  the  prob¬ 
lem  of  simulating  interstorm  runoff  dynamics  such  as  evapotranspiration  and 
long-term  base  flow.  Assuming  that  overland  flow  will  dominate  the  runoff 
process,  a  typical  kinematic  wave  model  will  be  used  to  model  a  watershed 
(Fig.  1).  The  watershed  has  two  overland  flow  planes  both  with  width  W, 
and  overland  flow  lengths  and  slopes  equal  to,  respectively,  li,  I2,  and  Sqi, 
5o2*  Rainfall  excess  is  assumed  to  be  uniformly  distributed  over  each  over¬ 
land  flow  plane.  The  physical  properties  and  the  initial  conditions  for  the 
watershed  represented  by  the  overland  flow  planes  are  assumed  to  be  uiu- 
f(nm  both  horizontally  and  vertically.  Rainfall  excess  is  calculated  by  first 


subtracting  an  initial  surface  loss 

r(r)  =  0  for  Pit)  ^  r  ^  0 . . (1) 

Kt)  =  ro(f)  for  P(r)  >  /„  f  2:  0 . (2) 

where  Pit)  =  the  cumulative  precipitation  over  the  watershed;  rit)  =  the 
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RG.  1.  Schematic  of  Watershed  Model 


rainfall  intensity  adjusted  for  surface  losses;  /  =  the  time  since  die  start  of 
rainfall;  ro(/)  =  the  measured  rainfall  intensi^;  and  =  the  depth  of  surface 
loss  assumed  to  be  uniform  over  the  watershed.  Second,  an  infiltration  loss 
is  calculated  by  the  Green  and  Ampt  (GA)  approach 


n  KXi 


i>K, 


(3) 


where  dF/dt  —  i(0  is  equal  to  the  infiltration  rate;  F  =  the  cumulative  in¬ 
filtration;  K,  »  the  soil’s  saturated  hydraulic  conductivity;  fl  »  the  product 
of  the  wetting  front  suction,  «(>/,  and  the  soil  volumetric  deficit  at  die  begin¬ 
ning  of  the  storm,  A0  =  ^  -  0*  <J>  =  the  soil’s  porosity;  and  0,  *=  the  initial 
volumetric  water  content.  The  GA  equation,  as  originally  developed,  is  <mly 
strictly  applicable  to  a  unifom  moisture  condition  at  the  soil  si^ace  or,  in 
the  case  of  rainfall  infiltxation,  a  ponded  surface  condition.  Modifications 
were  made  as  suggested  by  Mein  and  Larson  (1973)  and  Morel-Seytoux 
(1980)  to  use  the  GA  equation  for  unponded  surface  conditions  and  variable 
rainfall  rates. 

The  rainfall  excess  is  routed  overland  to  the  channel  using  the  kinematic 
wave  equations  of  motion 


dq  dy 

—  +  —  —  e  initial  condition:  r  =  0,  y  =  0 
dx  dt 


boundary  condition:  x  =  0,  y  =  0 . (4) 

9  =  0/' . (5) 


where  q  =  9(x,r)  is  the  flow  per  unit  width;  a  =  1.49So^/N;  So  =  the  slope 
along  die  overland  flow  length;  ^  =  a  roughness  coefficient  for  overland 
flow;  m  =  1.67;  y  =  y(x,r)  is  the  depth  of  flow;  and  x  -  the  distance  along 
the  direction  of  overland  fiow.  The  equations  will  be  solved  using  a  finite 
difference  scheme  proposed  by  Leclerc  and  Schaake  (1973). 

Flow  entering  the  channel  will  be  concentrated  immediately  to  the  wa¬ 
tershed  outlet  to  obtain  the  total  outflow,  Q.  Channel  routing  cf  the  flow 


using  the  kinematic  wave  equations  is  a  standard  procedure.  However,  the 
channels  of  the  test  catchments  used  in  this  research  are  short  and  have  little 
storage.  Performing  channel  routing  would  add  little  in  simulation  caq)ability, 
while  adding  greatly  to  the  computational  burdens  in  performing  the  uncer¬ 
tainty  analysis.  Consequently,  a  reasonable  assumption  was  made  that  chan¬ 
nel  routing  could  be  ignored. 

Parameter  Estimation 

Estimating  Green  and  Ampt  Parameters  from  Simple  Indices 
Rawls  et  al.  (1981,  1982,  1983)  and  Rawls  and  Brakensiek  (1982)  re¬ 
ported  on  the  relationship  between  simple  soil  indices  and  the  Green  and 
Ampt  parameters.  The  general  procedure  involved  has  been  to  relate  a  simple 
soil  index  such  as  texture  class  to  the  Brooks  and  Corey  parameters  (Rawls 
et  al.  1982),  porosity  4>*  residual  saturation,  6„  bubbling  pressure  i(f»,  and 
the  pore-size  distribution,  X.  These  parameters  are  used  to  describe  the  soU 
matric  suction  curve  by  the  relationship 


e  -  V 
<|)  -  e, 


(6) 


where  =  the  effective  saturation;  and  6  =  the  volumetric  water  content 
at  matric  suction  t|i.  The  Brooks  and  0)rey  parameters  are  then  used  to  cal¬ 
culate  the  wetting  front  sucdon,  by 

tly  =  <!»•  +  7^  ~  . (7) 

1  —  n 


where  =  the  water  entry  pressure;  *=  the  water  content  corresponding 
to  die  initial  soil  moisture  content  of  die  soil  prior  to  ponded  infUtradon; 
and  die  water  entry  pressure,  and  n  are  defined  as 


n  =  3X  +  2 . (9) 

Assuming  that  the  initial  moisture  content  is  equal  to  the  residual  saturation, 
die  formula  finally  derived  by  Brakensiek  (1977)  and  applied  by  Rawls  et 
al.  (1982)  is  obtained  as 

=  . . (10) 

#1-1 


Since  initial  moisture  conditions  will  be  estimated  using  an  API  index,  4*/ 
will  be  calculated  via  Eq.  7  in  model  simulations  of  infiltration.  The  satu- 
rrued  hydraulic  conductivity  is  calculated  using  Brutsaert’s  (1967)  solution 
of  die  (Zhilds  and  CoUis-George  permeability  integral 


4iJ(X  +  1)(X  +  1) 


(cm/s) 


(11) 


where  6,  =  —  6,.  The  parameter  a  was  set  equal  to  21.0  by  fitting  Eq. 
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11  to  an  average  of  observed  saturated  hydraulic  conductivity  values  for  a 
full  range  of  soil  texture  classes  (Rawls  et  al.  1982).  Rawls  and  Brakensiek 
(1989)  developed  the  regression  relationships  shown  here  between  simple 
soil  indices;  p^icle  size  distribution  and  total  porosity,  and  the  Brooks  and 
Corey  parameters. 

The  regression  equations  (Rawls  and  Brakensiek  1989)  are 
e,.  =  -0.0182482  +  0.00087269(5A)  +  0.005 13488(CL)  +  0  02939286(P/?) 

-  0.00015395(CL^)  -  0.0010827(5A)(P/?)  -  0.0001 8233(CL^)(PP^) 

+  0.00030703(CL^)(P/?)  -  0.0023584(P/f^)(CZ-) . (12) 

X  =  exp  [-0.7842831  +  0.0177544(SA)  -  1.062498(PP)  -  0.00005304(5A^) 

-  0.00273493(CL*)  +  1.11134946(PP^)  -  0.03088295(5A)(PP) 

+  0.00026587(SA^)(P/?^)  -  0.0061052(CL^)(PP^)  -  0.00000235(SA*)(CL) 

+  0.00798746(CL*)(PP)  -  0.0067449 1(PP*)(CL)] . (13) 

!  «=  exp  [5.3396738  +  0.1845038(CL)  -  2.48394546(PP) 

j  -  0.00213853(CL^)  -  0.04356349(SA)(PP)  -  0.61745089(CL)(PP) 

j  +  0.00143598(SA*)(PP^)  -  0.00855375(CL^)(PP^)  -  0.00001282(SA*)(CL) 

t+  0.00895359(CL*)(PP)  -  0.00072472(SA*)(PP)  +  0.0000054(CL*)(SA) 

+  0.50028060(P/?^)(CL)] . (14) 

where  PR  =  <|)  =  total  porosity  (vol/vol);  CL  =  %  clay  (particle  size);  SA 
=  %  sand  (particle  size);  ^i,  =  bubbling  pressure  (cm);  X  =  pore>size  dis¬ 
tribution;  and  6,  residual  saturation  (vol/vol).  The  data  for  these  regres- 
sitms  were  obtained  from  the  numerous  soil  surveys  performed  in  die  United 
States. 

The  data  needed  to  apply  these  regressions  can  be  conveniendy  obtained 
from  the  SOILS-5  data  base  compiled  by  the  Soil  Conservadon  Service  and 
maintained  by  the  Corps  of  Engineers  (SOILS-S  1983).  The  relevant  infor- 
madon  provided  by  these  data  is  the  pardcle-size  distribudon,  the  percent 
sand,  clay,  and  silt,  and  the  porosity  for  different  horizons  of  a  soil  series. 

Selection  of  Soil  Parameter  Probability  Distributions 
In  the  case  of  the  watershed  model  used  in  this  research,  the  chosen  pa¬ 
rameters*  probability  distribudon  should  describe  the  chance  that  the  param¬ 
eter  will  have  a  certain  value  for  a  given  storm  event.  The  random  variadon 
of  a  particular  infiltration  parameter  between  events  might  be  affected  by 
biolq^c  (macropores),  cultural  (management  practice),  seasonal  (crop  giowtii), 
and  meteorologic  (rain-induced  surface  crusting)  factors. 

The  distribution  of  soil  properties  within  a  given  texture  class  or  for  a 
given  particle-size  distribution  does  not  give  information  on  the  random  vari- 
'  ation  of  soil  hydraulic  properties  in  the  field.  Rather,  the  distribution  of  in¬ 
filtration  parameter  values  for  a  given  soil  property  index  is  a  measure  of 
the  uncertainty  in  predicting  the  infiltration  properties  of  a  field  sample  using 
that  index.  The  correspondence  between  the  probability  distribution  derived 
from  a  simple  index  and  a  field  measurement  program  would  be  coincidental 
at  best. 
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TABLE  1.  ValtiM  for  Watershed  Model's  Stochastic  Parameters 


Watershed 

(1) 

Sand 

Clay 

Porosity 

Mean 

(%) 

(2) 

Range 

(%) 

(3) 

Mean 

(%) 

(4) 

Range 

{%) 

(5) 

Mean 

(vol/vol) 

(6) 

Coefficient 
of  variation 
(7) 

8 

3.4 

0.4-6.4 

23.9 

19.5-28.4 

■3B 

0.1 

9 

3.1 

0.3-5.9 

21.4 

'16.6-26.2 

0.1 

10 

2.5 

0.0-5.0 

21.7 

16.9-26.5 

0.51 

0.1 

11 

2.5 

0.0-5.0 

23.1 

18.5-27.7 

0.51 

0.1 

12 

2.5 

0.0-5.0 

24.6 

20.2-29.0 

0.51 

0.1 

19 

2.5 

0.0-5 .0 

26.4 

22.2-30.6 

0.51 

0.1 

20 

2.5 

0.0-5.0 

27.4 

23.4-31.4 

0.51 

0.1 

The  probability  distribution  obtained  from  a  simple  index  might  be  useful 
if  the  ^stribution  of  the  parameter's  values  (e.g.,  the  values  of  sahirated 
hydraulic  conductivity)  bounds  the  range  of  the  values  for  the  actual  field 
parameter.  Furthermore,  the  parameter  distribution  obtained  may  be  even 
more  useful  from  an  engineering  perspective  if  it  not  only  bounded  the  po> 
tential  range  of  the  field  parameter,  but  also  led  to  a  conservative  watershed 
model  prediction  of  runoff.  For  example,  the  simple  index  probability  dis¬ 
tribution  should  overpredict  the  chance  that  the  saturated  hydraulic  conduc¬ 
tivity  will  be  less  than  the  true  field  value,  which  in  turn  will  result  in  an 
overprediction  of  the  chance  that  runoff  will  exceed  a  certain  value.  The 
overprediction  resulting  from  the  simple  index  may  be  considered  to  provide 
a  safety  factor  that  would  be  useful  in  design  situations. 

The  probability  distribution  for  the  Green  and  Ampt  infiltration  parameters 
will  be  derived  from  assumed  distributions  for  particle  size  distribution  and 
total  porosity  via  Eqs.  7,  11,  12,  13,  and  14.  The  SOILS-5  data  base  pro¬ 
vides  die  range  in  particle-size  distribution  and  total  porosity  for  the  wa¬ 
tershed  soils  (see  Table  1).  The  particle-size  distribution,  as  detined  by  the 
percentages  of  sand  and  clay,  was  assumed  to  be  distributed  uniformly  over 
this  range.  The  total  porosity  obtained  from  the  SOlLS-5  data  base  was  as¬ 
sumed  to  have  a  nonnal  distribution  with  a  coefficient  of  variation  equal  to 
0.1  in  analogy  to  that  found  for  variation  within  a  texture  class  (McCuen  et 
al.  1981). 

Deriving  the  distributions  of  the  Green  and  Ampt  parameters  K,  and  i|f/  is 
an  intermediate  step  in  deriving  the  volume  and  peak  discharge  probability 
distributions.  The  probability  distribution  for  K,  is  independent  of  antecedent 
moisture  and  was  derived  using  Monte  Carlo  simulation.  The  mean,  standard 
deviation,  skew,  and  kurtosis  derived  by  Monte  Carlo  simulation  ate  shown 
in  Table  2. 

The  distribution  of  is  conditional  on  the  antecedent  moisture  condition, 
and  is  not  easily  repotted  in  general.  However,  the  i|>/  corresponding  to  tiie 
mean  infiltration  parameters  and  initial  moisture  content  equal  to  the  residual 
saturation  is  shown  in  Table  3  for  informational  purposes.  The  derived  mean 
values  for  K,  and  <}i/  will  not  be  equivalent  to  those  calculated  directly  from 
average  values  of  total  porosity  and  particle  distribution  as  can  be  seen,  for 
example,  by  comparing  the  values  of  shown  in  Tables  2  and  3. 

No  attempt  will  be  made  to  assess  the  sensitivity  of  model  predictions  to 
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TABLE  2.  Derived  Values  for  Saturated  Hydraulic  Conducttvtty 


Watershed 

(1) 

Mean 

(in./hr) 

(2) 

Standard  deviation 
(in./hr) 

(3) 

Coefficient 
of  variation 
(4) 

Skew 

(5) 

Kurtosis 

(6) 

8 

MSM 

0.044 

MEM 

1.6 

6.6 

9 

■EH 

0.045 

1.6 

6.5 

10 

0.058 

0.044 

1.6 

6.2 

11 

0.055 

0.042 

1.6 

6.4 

12 

0.051 

0.041 

0.80 

1.7 

6.6 

19 

0.047 

0.039 

0.83 

1.7 

7.0 

20 

0.045 

0.038 

0.84 

1.8 

7.2 

Note:  1.0  in.  =  2.54  cm,  1.0  in./hr  =  2.54  cm/h. 


the  assumed  probability  distributions  (i.e.,  normal  or  uniform)  for  particle- 
size  distribution  or  total  porosity.  In  verifying  the  model,  the  focus  will  be 
(m  con^)aring  model  predictions  within  one  and  two  standard  deviations  about 
the  mean  prediction.  Presumably,  the  standard  deviation  about  the  mean  pre¬ 
diction  will  be  most  affected  by  the  standard  deviation  assumed  for  model 
parametos.  The  tails  of  the  assumed  distributions,  and  consequently  the  shape, 
should  be  less  important.  Certainly,  the  best  way  to  test  sensitivity  of  model 
pr^ctions  to  assumed  parameter  distributions  is  to  try  some  different  dis¬ 
tributions.  However,  this  was  not  possible  given  the  computational  burdens 
required  to  derive  model  predictions  by  Monte  Carlo  simulation. 

Estimaticn  of  Initial  Abstraction 

The  initial  abstraction  is  defined  as  the  water  that  is  not  free  to  flow  over¬ 
land  to  a  stream.  These  losses  are  due  to  the  interception  of  rain  by  crops 
and  the  depression  storage  due  to  the  microrelief  of  the  surface.  Viessman 
et  al.  (1977)  summarize  studies  that  could  be  used  to  estimate  interception 
loss.  However,  estimation  of  the  depression  storage  from  soil  survey  infor¬ 
mation  is  extremely  difficult.  Linden  (1979)  estimated  depression  storage  for 
cultivated  soils  using  a  microrelief  model.  The  key  parameters  in  this  model 
are  random  roughness  and  land  surface  slope  (random  roughness  is  essen¬ 
tially  a  measure  of  the  variation  of  soil  heights  from  the  plane  that  would 


TABLE  3.  Deterministic  Vaiues  for  Watershed  Green  and  Ampt  and  Brooks  and 
Corey  Parameters 


Watershed 

(1) 

Sand 

(%) 

(2) 

Clay 

(%) 

(3) 

Porosity 

4 

(voi/voQ 

(4) 

Residual 

saturation, 

•r 

(ven/vol) 

(5) 

ERsctiva 

saturation 

e, 

(vol/voO 

<6) 

Wetting 

front 

suction 

"lyPn.) 

(7) 

Saturated 

hydraulic 

conductivity 

X,(in./hr) 

(8) 

BubMng 

pressure 

4>»(in.) 

(9) 

Pore  size 
distitxitloo 

X 

(10) 

t 

3.4 

23.9 

■51 

0.080 

0.430 

18.62 

0.0430 

24.19 

0.296 

9 

3.1 

21.4 

ESI 

0.074 

0.436 

18.14 

0.0487 

23.67 

0.306 

10 

2.5 

21.7 

0.075 

0.435 

18.41 

0.0468 

24.03 

0.305 

11 

2.5 

23.1 

OJl 

0.07S 

0.432 

18.69 

0.0435 

24.35 

0.299 

12 

2.5 

24.6 

0.51 

0.082 

0.428 

19.02 

0.0400 

24.72 

0.293 

19 

2.5 

26.4 

0.51 

0.086 

0.424 

19.45 

0.0364 

25.19 

0.286 

20 

2.5 

27.4 

0.51 

0.068 

0.422 

19.69 

0.0342 

25.48 

0.281 

Mole:  1.0  in.  «  2.S4  cm,  1.0  fai./hr  •  2.S4  cm/b. 
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define  zero  cut  and  fill).  Application  of  Linden’s  results  are  hampered  be¬ 
cause,  as  Allmaras  et  al.  (19^)  found,  random  roughness  in  an  agricultural 
soil  is  a  function  of  the  plowing  tool,  initial  water  content  prior  to  plowing, 
and  soil  preparation.  Consequently,  random  roughness  and,  correspondingly, 
depression  storage  are  very  difficult  to  estimate  from  soil  survey  informa¬ 
tion. 

Linden’s  results  indicated  that  for  the  surface  slopes  of  the  test  watersheds, 
3-7%,  the  depression  storage  could  vary  between  0.08  and  0.40  in.  (0.20 
and  1.05  cm),  corresponding  to  random  roughness  values  between  0.3  and 
1.6  in.  (0.8  and  4.0  cm)  (Linden  1979).  Combining  this  depression  storage 
with  potential  interception  loss  would  lead  to  a  surface  loss  that  at  the  high 
end  could  account  for  the  observed  losses  in  80%  of  the  observed  events. 

Given  that  there  are  no  simple  indices  for  estimating  random  roughness 
from  soil  survey  information  and  the  large  impact  that  assumptions  about 
this  parameter  can  have  on  predicted  runoff,  the  assumption  was  made  that 
the  initial  loss  is  zero.  This  assumption  should  lead  to  a  model  bias  towards 
overprediction  on  the  average.  The  actual  bias  can  be  observed  by  comparing 
mean  predicted  runoff  with  observed  runoff  values.  Whether  the  bias  is  in 
a  range  reasonably  attributed  to  surface  loss  can  be  judged  in  the  context  of 
future  research.  The  advantage  of  this  assumption  is  that  the  surface  loss 
does  not  end  up  being  a  “fudge  factor”  that  is  varied  to  account  for  model 
bias. 

Kinematic  Wave  Parameter  Estimates 

The  estimates  of  overland  flow  lengths  and  slopes  for  overland  flow  planes 
were  obtained  from  2-ft  contour  interval  topographic  maps.  The  roughness 
values  for  overland  flow  were  estimated  as  a  best  average  value  for  all  runoff 
events  based  on  published  values  [summary  given  by  Hjemfelt  (1986)).  The 
values  for  these  parameters  are  given  in  Table  4  for  each  test  watershed. 

Estimating  Initial  Moisture  Content 

An  estimate  of  the  average  moisture  conditions  prior  to  a  rainfall  event  is 
needed  to  implement  the  Green  and  Ampt  method.  The  moisture  condition 
was  estimated  by  developing  regression  relationships  between  an  antecedent 
precipitation  index  (API)  and  essentially  point  measurements  of  soil  mois¬ 
ture.  The  regression  relationships  will  then  be  used  to  predict  the  average 
moisture  condition  prior  to  a  rainfall  event  given  the  API. 

There  are  a  number  of  difficulties  with  this  approach.  First,  the  use  of 
soil  moisture  data  to  develop  this  relationship  is  a  departure  from  the  un¬ 
gauged  analysis  procedure.  However,  soil  moisture  measurements  are  rarely 
available,  and  it  was  felt  that  investigating  their  usefulness  was  an  important 
opportunity.  Second,  the  assumption  is  made  that  the  estimates  of  soil  mois¬ 
ture  obtained  from  grab  samples  (point  measurements)  are  indicative  of  the 
general  watershed  moisture  condition.  Certainly,  there  will  be  spatial  vari¬ 
ation  of  soil  moisture  across  a  watershed.  However,  the  point  soil  moisture 
is  the  best  indicator  of  the  average  watershed  moisture  condition  available, 
although  it  has  drawbacks. 

Third  and  finally,  the  assumption  that  an  API  can  be  used  to  assess  av¬ 
erage  watershed  conditions  ignores  the  effect  of  evapotranspiration  between 
storm  events  on  available  moisture.  Again,  API  is  not  an  ideal  proxy  for 
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TABLE  4.  Kinematic  Wave  Parameters  for  Test  Watersheds 


Watershed 

(1) 

Area* 

(acre) 

(2) 

Overland  flow 
length  /  (ft) 

(3) 

Overland 
slope  So 
(4) 

Roughness 

coefficient 

N 

(5) 

initial 

abstraction 
/.  (in.) 

(6) 

8 

3.84 

280 

0.07 

0.1 

0 

8 

3.64 

220 

0.07 

O.I 

0 

9 

3.93 

280 

0.04 

0.1 

0 

9 

4.02 

220 

C  05 

0.1 

0 

10 

4.16 

250 

0.03 

0.1 

0 

10 

4.01 

200 

0.05 

0.1 

0 

11 

4.16 

210 

0.04 

0.1 

0 

11 

4.26 

200 

0.03 

0.1 

0 

12 

3.93 

280 

0  03 

0.1 

0 

12 

3.97 

270 

0.02 

0.1 

0 

19 

3.62 

390 

0.07 

0.1 

0 

19 

3.62 

280 

0.04 

0.1 

0 

20 

2.48 

500 

0.05 

0.1 

0 

20 

2.48 

220 

0.05 

0.1 

0 

'Area  after  19S9. 

Note:  1.0  acre  =  0.40S  'la;  1.0  ft  =  3.28  m,  1.0  in.  =  2.54  cm. 


soil  moisture;  l>ut  for  ungauged  analysis  with  a  single-event  model,  it  is 
probably  the  only  one  available. 

The  soil  moisture  data  were  obtained  at  various  unrepoited  locations  within 
the  test  watersheds  {The  Centred  no  date).  The  data  were  obtained  by  making 
gravimetric  measurements  of  soil  moisture  for  specimens  taken  ^m  fist¬ 
sized  field  samples.  The  sampling  interval  was  unevenly  spaced  over  about 
30  years  of  the  measurement  program.  The  volumetric  water  content  was 
calculated  from  die  reported  gravimetric  values,  6^  (fraction  of  dry  unit  weight), 
by  applying  the  following  formula: 

0  =  =  1.362200 . (15) 

Iw 

where  0  =  the  volumetric  water  content;  yn,  =  the  unit  weight  of  water;  and 
750  =  the  dry  unit  weight  of  soil,  taken  as  85  Ib/cu  ft  for  all  soils  in  the 
study  area. 

The  soil  measurements  are  most  plentiful  for  the  months  during  the  grow¬ 
ing  season,  April-October.  The  measurements  were  taken  at  depths  between 
0.0  and  6.0  ft.  For  the  purposes  of  calculating  infiltration,  soil  moisture 
measurements  from  0.0  to  1.0  ft  were  considered  most  important.  These  soil 
moisture  measurements  covered  a  full  range  of  soil  moisture  conditions.  A 
linear  regression  relationship  was  established  between  the  5-,  15-,  and  30- 
day  API  and  the  estimated  soil  moisture  to  explore  the  effects  of  short-, 
m^um-,  and  long-term  antecedent  rainfall.  Both  multiple  icgression  and 
simple  regressions  were  explored,  as  can  be  seen  from  the  results  shown  in 
Fig.  2.  The  multiple  coefficient  of  determination,  is  not  tremendously 
impressive  for  any  of  the  months  examined.  The  results  indicate  that  the 
best  simple  regression  was  not  significantly  poorer  than  the  multiple  regies- 


FIG.  2.  Multiple  Linear  Regression  between  5-,  15-,  and  30-Day  API  and  Soil 
Moisture  at  0-1  ft  Depth 


sion  relationships.  Consequently,  the  best  simple  regression  relationship  was 
used  for  predicting  average  watershed  moisture  conditions. 

The  error  made  in  estimating  the  average  moisture  condition  via  API  is 
computed  by  using  the  final  ^egression  results  shown  in  Table  5.  The  esti¬ 
mated  error  is  modeled  by  a  normal  distribution  whose  mean  value  is  given 
by  the  regression  for  each  month  and  standard  deviation  is  estimated  by  the 
standard  error  of  the  regression.  Technically,  the  error  distribution  should 
be  proportional  to  the  product  of  the  standard  error  and  the  student’s  /-dis¬ 
tribution.  However,  the  number  of  degrees  of  freedom  available  for  the 
regression  makes  the  difference  between  the  normal  and  /-distributiDn  small 
enough  to  be  ignored. 


TABLE  5.  Simple  Regression  Equations  for  API  and  Soil  Moisture  Using  Best 
Index 


Month 

(1) 

Abscissa 

intercept 

a 

(2) 

Slope 

b 

(3) 

Simple 

correlation, 

R 

(4) 

Number  of 
API  days* 
(5) 

Mean 

predicted 

moisture 

content 

(vol/vol) 

(6) 

Standard 

error 

(vol/vol) 

(7) 

Number  of 
observations 
(8) 

4 

0.32 

0.022 

0.33 

15 

0.35 

0.07 

230 

S 

0.34 

0.004 

0.04 

5 

0.34 

0.07 

189 

6 

0.26 

0.046 

0.71 

15 

0.40 

0.11 

227 

7 

0.14 

ESHl 

0.58 

30 

0.21 

0.07 

171 

8 

0.23 

0.06< 

0.38 

15 

0.33 

199 

9 

0.16 

0.043 

0.77 

30 

0.32 

WSSM 

212 

10 

0.23 

0.036 

0.47 

30 

0.30 

0.10 

188 

*6<  ~  a  4-  biP),  where  /*  is  the  antecedent  precipitation  for  days  designated  in  column  5. 
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Stochastic  Analysis  Results 

Stochastic  model  predictions  were  made  by  using  Monte  Carlo  simulation 
to  dmve  the  probability  of  runoff  volumes,  peak  discharge,  and  time-to- 
peak  discharge  for  a  given  rainfall  event.  The  distributions  were  estimated 
by  performing  sufficient  Monte  Carlo  simulations  to  obtain  stable  estimates 
of  the  distribution’s  mean,  standard  deviation,  skew,  and  kurtosis.  The  dis¬ 
tributions  were  estimated  for  101  events  observed  for  tliv  test  watersheds. 
Other  rainfall  events  were  available  from  the  test  watershed  records,  but 
various  observation  errors  excluded  their  use. 

The  statistics  for  the  predicted  volume  and  peak  discharge  stabilized  in  a 
reasonable  number  of  simulations  (5,000).  However,  the  derived  standard 
deviation,  skew,  and  kurtosis  for  the  time-to-peak  discharge  did  not  stabilize 
as  readily,  requiring  at  least  50,000  simulations.  Since  the  estimation  of 
tinie-to-peak  disdiatge  statistics  led  to  onerous  computer  time,  only  die  mean 
time-to-peak  discharge  was  calculated. 

The  model  was  verified  by  comparing  both  model  mean  runoff  predictions 
and  predicted  exceedance  probabilities  as  determined  from  derived  distri- 
tnitions  with  corresponding  runoff  observations  and  exceedance  probabilities 
estimated  from  observations.  A  comparison  of  mean  predicted  and  observed 
runoff  resulted  in  a  bias  towards  overprediction  for  runoff  volumes  of  0.20 
in.  and  for  peak  discharge  ah  overprediction  of  about  1.3  cfs.  The  prediction 
.of  time-to-peak  discharge  waas  biased  towards  an  underprediction  of  about 
24  minutes  on  the  avoage.  However,  any  individual  prediction  of  volume, 
peak  discharge,  or  time-to-peak  discharge  can  differ  from  the  observed  value 
by  a  much  larger  amount  ^an  the  average  bias,  as  can  be  seen  from  Hgs. 
3-5. 

The  exceedance  probability  comparisons  were  made  at  plus  or  nunus  one 
or  two  standard  deviations  for  the  derived  distributions  (see  Fig.  6).  The 
observed  exceedance  or  nonexceedance  probabilities  are  defined  as,  respec- 


FICL  3.  Predicted  Mean  Volumes  versus  Observed  Volumes 
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FIG.  4.  Predicted  Mean  Peak  Discharge  versus  Observed  Volumes 


tively,  fbe  firaction  of  observed  occurrences  that  are  greater  than  the  mean 
prediction  plus  one  or  two  standard  deviations,  and/or  less  than  the  mean 
prediction  minus  one  or  two  standard  deviations. 

To  estimate  die  predicted  ^ceedance  probability,  a  functional  form  for 
die  derived  distributions  was  determined  by  comparing  the  relationship  of 
the  coeificient  of  variation  of  the  skew  or  kurtosis  for  the  delved  runoff 
distributions  with  that  of  the  normal  and  lognormal  distributions,  as  shown 


FKau  6.  Prsdlded  Mean  Time  to  Peak  Discharge  versus  Observed  Peak  DIs* 
charge 
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POF 


-Meon  Volue 


Watershed  Model  Prediction  Variable 


FIG.  6.  Definition  of  Error  Bounds  and  Corresponding  Error  Exceedance  Prob- 
ablinies 


i&  Hgs.  7  and  8.  A  comparison  of  the  simulation  results  shows  that  the 
derived  distributions  for  both  runoff  volume  and  peak  discharge  are  better 
represented  by  a  normal  distribution  than  a  lognormal  distribution. 

Assuming  ^at  the  derived  distributions  ate  normal,  model-predicted  ex¬ 
ceedance  probabilities  at  plus  or  minus  one  and  two  standard  deviations  about 
the  mean  are  equal  to  16%  and  2.5%,  respectively.  The  estimated  observed 
exceedance  probabilities  for  the  volume  and  peak  discharge  analysis  are 


RQ.  7.  Comparison  of  Derived  Skew  for  Runoff  Voiumes  and  Peak  Discharges 
with  Those  for  Normal  and  Two-Parameter  Lognormal  Distributions 
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FIG.  8.  Comparison  of  Derived  Kurtosis  for  Runoff  Volumes  and  Peak  Dis¬ 
charges  with  Those  for  Normal  and  Two-Parameter  Lognormal  Distributions 

compared  to  the  predicted  ones  in  Table  6.  The  comparisons  indicate  that 
die  model  conservatively  predicts  exceedance  probability  at  plus  one  stan¬ 
dard  deviation  for  both  prediction  variables  and  at  two  standard  deviations 
for  runoff  volume.  However,  the  model  prediction  of  discharge  is  not  con¬ 
servative  at  two  standard  deviations  in  that  the  estimated  observed  exceed¬ 
ance  probability  is  3.0%,  whereas  the  predicted  one  is  2.5%.  Consequently, 
the  results  demonstrate  that  the  model  is  conservative  at  one  standard  d«^ 
viation,  but  not  necessarily  at  two  standard  deviations.  A  corresponding 
comparison  of  nonexceedance  probabilities  also  shown  in  Table  6  indicates 


TABLE  6.  Comparison  of  Observed  and  Predicted  Exceedance  Probabilities 
(Based  on  101  Events)   


Number  of 
standard 
deviations 
(1) 

Exceedance  (%)• 

Nonexceedance  (%)*’ 

Predicted 

(2) 

Observed 

(3) 

Predicted 

(4) 

Obsenred 

(5) 

(a)  Volume 

1.0 

16.0 

l.O 

16.0 

44.6 

2.0 

2.5 

2.5 

12.9 

(b)  Peak  Discharge 

1.0 

16.0 

9.9 

16.0 

31.7 

2.0 

2.5 

3.0 

2.5 

9.9 

TExceedance  probability  «=  the  chance  that  ninoff  volume  or  peak  discharge  will  exceed 
die  mean  prediction  by  the  given  number  of  standard  deviations. 

^onexceedance  probability  «  the  chance  that  runoff  volume  or  peak  discharge  will  be 
less  than  die  mean  minus  the  given  number  of  standard  deviations. 
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that  the  model  underpredicts  the  chance  that  the  runoff  will  be  less  than  one 
or  two  standard  deviations.  Thus,  the  model  is  conservative  in  this  respect. 

Conclusions 

The  uncertainty  in  model  predictions  was  conservatively  estimated,  for  the 
most  part,  by  considering  ^e  error  resulting  from  using  simple  indices  to 
estimate  infiltration  parameters.  In  other  words,  the  assumed  error  model  for 
inflltration  parameter  estimates  resulted  in  estimated  model  prediction  errors 
that  were  somewhat  greater  than  the  observed  prediction  error.  This  con¬ 
servatism  may  be  useful  in  that  it  could  be  used  to  obtain  a  safety  factor  in 
a  design  situation. 

The  conservatism  of  the  estimated  prediction  uncertainty  might  be  attrib¬ 
uted  to  the  assumptions  made  in  choosing  parameter  probability  distributions 
and  in  assuming  zero  surface  loss.  The  conservatisin  in  the  estimated  pre¬ 
diction  uncertainty  would  be  lessened  if  some  type  of  peaked  distribution 
were  used  instead  of  the  uniform  distribution  assumed  for  particle-size  dis¬ 
tribution.  Estimation  of  a  surface  loss  may  have  helped  in  reducing  the  bias 
towards  overprediction  of  runoff.  However,  further  research  in  quantifying 
surface  losses  for  ungauged  watersheds  would  be  helpful  in  attempts  to  es¬ 
timate  this  parameter  for  the  model.  Additionally,  the  infiltration  estimation 
procedure  could  be  improved  by  finding  a  better  indicator  of  soil  moisture 
'  than  an  antecedent  precipitation  index. 

The  estimation  procedure  used  for  infiltration  parameters  could  be  useful 
when  ^>plying  a  PBS  type  model  to  derive  flow-fiequency  curves  (Eagleson 
1972).  However,  application  of  this  estimation  prcMxdure  to  more  gauged 
basins  would  be  necessary  to  determine  whether  it  is  consistently  conser¬ 
vative  and  to  determine  the  degree  of  conservatism. 
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Appendix  II.  Notation 

The  following  symbols  are  used  in  this  paper: 

F  =  cumulative  infiltration  (in.); 
la  =  initial  abstraction  (in.); 

I  =  infiltrat: 'a  rate  (in./hr); 

Kf  =  saturated  hydraulic  conductivity  (in. /hr); 

/  =  length  along  overland  flow  plane  (ft); 

P  =  cumulative  precipitation  (in,); 

Q  =  total,  watershed  outflow  (cfs); 
q  —  flow  per  unit  width  of  channel  (cfs/ft); 
m  =  kinematic  wave  discharge  exponent  =  1.67; 

N  =  roughness  coefficient  for  overland  flow; 
n  =  3X  +  2; 

r  =  rainfall  intensity  adjusted  for  surface  loss  (in./b*'}; 
To  =  rainfall  intensity  (in. /hr); 

So  =  slope  along  ove;.:  id  ^ow  plane  (ft/ft); 

St  =  effective  saturation; 
t  =  time  (sec); 
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W  ~  width  of  overland  flow  plane  (ft); 

X  -  distance  along  overland  flow  plane  (ft); 
y  =  flow  drath  along  overland  flow  plane  (ft); 
o  =  \A9S]r/N\ 

"Isd  -  soil  dry  unit  weight  (Ib/cu  ft); 
yw  =  unit  weight  of  water  (Ib/cu  ft); 

A6  =  volumetric  water  deficit  (vol/vol); 

0  =  voluntetric  water  content  (vol/vol); 

6o  =  gravimetric  water  content  (fraction  of  dry  unit  weight); 
0,  =  ^  -  0r  (vol/vol); 

0,  =  initial  water  content  (vol/vol); 

0r  =  residual  saturation  (vol/vol); 

X  =  pore'S.ze  distribution; 
a  =  standard  deviation; 

<t>  =  porosity  (vol/vol); 

^  =  matric  suction  (in.); 

=  bubbling  pressure  (in.); 

=  wetting  ftont  suction  (in.); 

»|»w  =  water  entry  pressure  (in.);  and 
ft  =  0®-)- 

i^PENDix  III.  Conversion  to  SI  Units 


To  convert 

To 

Multiply  by 

acre 

ha 

0.405 

ft 

cm 

30.48 

in. 

cm 

2.54 

Ib/cu  ft 

N/m’ 

157.05 

cu  ft/sec 

m’/s 

0.0283 
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